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f\j Abstract 



Polyhedra form an established abstract domain for inferring runtime 

properties of programs using abstract interpretation. Computations 

on them need to be certified for the whole static analysis results to 

CO be trusted. In this work, we look at how far we can get down the 

^^ road of a posteriori verification to lower the overhead of certification 

1 of the abstract domain of polyhedra. We demonstrate methods for 

n , making the cost of inclusion certificate generation negligible. From a 

• performance point of view, our single-representation, constraints-based 

^ implementation compares with state-of-the-art implementations. 

^_^ In static analysis by abstract interpretation [:6||, sets of reachable states, 

J> which are in general infinite or at least very large and not amenable to 

■^ tractable computation, are over-approximated by elements of an abstract 

QQ domain on which the analyzer applies forward (resp. backward) steps 

(^ corresponding to program operations (assignments, tests. . . ) as well as 

■^ "joins" corresponding to control points with several incoming (resp. outgoing) 

^^ edges. When dealing with numerical variables in the analyzed programs, one 

, of the simplest abstract domains consists in keeping one interval per variable, 

:• and the forward analysis is known as interval arithmetic. Interval arithmetic 

. ^H however does not keep track of relationships between variables. The domain 

^ of convex polyhedra \7\ tracks relationships of the form ^^ ajXj cxi b where 

^ the tti and b are integer (or rational) constants, the Xi are rational program 

variables, and M is <, < or =. 

The implementor of an abstract domain faces two hurdles: the implemen- 
tation should be reasonably efficient and scalable; it should be reasonably 
bug-free. As an example, the Parma Polyhedra Library (PPL) f\], version 1.0, 
which implements several relational numerical abstract domains, comprises 
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260,000 lines of C++; despite the care put in its development, it is probable 
that bugs have slipped through. The same applies to the Apron library nU. 

Such hurdles are especially severe when the analysis is applied to large- 
scale critical programs (e.g. in the AsTREE system [5], targeting avionics 
software). For such systems, normal compilers may not be trusted, result- 
ing in expensive post-compilation checking procedures, and prompting the 
development of the CompCert certified compiler il2i: this compiler is pro- 
grammed, specified and proved correct in COQ [2l]. We wish to extend this 
approach to obtain a trusted static analyzer; this article focuses on obtaining 
a trusted library for convex polyhedra, similar in features to the polyhedra 
libraries at the core of PPL and Apron. 

One method for certifying the results of a static analysis is to store 
the invariants obtained by an untrusted analyzer at (roughly) all program 
points, then check that they are inductive using a trusted checker: each 
statement is then a Hoare triple that must be checked. Unfortunately, storing 
invariants everywhere proved impractical in the Astree analyzer due to 
memory consumption; we then opted to recompute them. Our (future) 
analyzer will thus store invariants only at loop heads, and thus, for control 
programs consisting of one huge control loop plus small, unrolled, inner loops, 
will store only a single invariant. It will then enter a checking phase which 
will recompute, in a trusted fashion, all intermediate invariants. Efficiency is 
thus important. 

The main contribution of our article is an efficient way of implementing 
a provably correct abstract domain of polyhedra. Efficiency is two-fold: 

1. In proof effort: most of the implementation consists in an untrusted 
oracle providing certificates of the correctness of its computations; 
only a much smaller certificate checker, consisting in simple algorithms 
(multiplying and adding vectors, replacing a variable by an expression), 
needs to be proven correct in the proof assistant. 

2. In execution time: the expensive parts of the computations (e.g. linear 
programming) are inside the untrusted oracle and may use efficient pro- 
gramming techniques unavailable in parts that need formal proofs. We 
do not compute certificates as an afterthought of polyhedral computa- 
tions: close examination of the algorithms implementing the polyhedral 
operators revealed that they directly expose the elements needed to 
build certificates. Simple bookkeeping alleviates the need to rebuild 
them after the fact. The overhead of making the operators certifying 
is thus very limited. This contrasts with earlier approaches [4] based 
on a posteriori generation of witnesses, which had to be recomputed 
from scratch using linear programming. 

A second contribution is a complete implementation of the abstract 
domain of polyhedra in a purely constraints-based representation. Most 



libraries used in static analysis, including PPL and Apron, use a double 
description: a polyhedron is both an intersection of half-spaces (constraints) 
or the convex hull of vertices, half-lines and lines (generators), with frequent 
conversions. Unfortunately, the generator representation is exponential in 
the number of constraints, including for common cases such as hypercubes 
(e.g. specification of ranges of inputs for the program). We instead chose 
to represent the polyhedra solely as lists of constraints, with pruning of 
redundant ones. Our implementation uses sparse matrices of rational numbers 
and uses efficient techniques for convex hull [20| and emptiness testing by 
linear programming [9]. 

We applied our library to examples of polyhedral computations obtained 
by running the Pagai static analyzer [lOj on benchmark programs. Despite 
a common claim that implementations based on the double representation 
are more efficient than those based on constraints only, our library reaches 
performance comparable to the Apron library together with the high-level 
of trust brought by our COQ certificate checker. 

The remainder of this paper is organized as follows. After having stated 



the conventions we are using (§1), we define correctness criteria for the 



operators of the abstract domain (§2), which all reduce to inclusion properties 
for which certificates are presented as Farkas coefficients (also known as 
Lagrange multipliers). Such certificates may also be cheaply generated for 



the convex hull (§4). Both forward step and convex hull operations reduce 



internally to a form of projection. Some design choices of our implementation 



are then described (§5), including how to keep the representation size of 



the polyhedra reasonable. Last before conclusion (§7), an experimental 



evaluation and accompanying results are presented (§6). 

1 Definitions and Notations 

In the remainder of this article, we use the following notations and definitions. 

1. C: a linear constraint of the form a ■ x < b where a is a vector of 
rational constants, 6 is a rational and x G Q" is the vector of the 
analyzed program variables. Such a linear constraint, or constraint for 
short, can be viewed as a half-space in an n-dimensional space. We 
write C = a ■ X > b for the complementary half-space. 

2. P: a convex polyhedron, not necessarily closed, represented as a set of 
constraints. We call "size of the representation of P" the number of 
constraints that P is made of. 

3. satisfaction: saying that point y of Q" satisfies a constraint C = a-x < b 
means that a ■ y < b. By extension, a point y satisfies (or is in) 
polyhedron P if it satisfies all of its constraints. We write this: Sat P y. 



Given that a constraint C can be regarded as polyhedron with only 
one constraint, we also write: Sat C y. 

4. Given our focus on the abstract domain of polyhedra we shall adopt 
the following vocabulary. 

(a) The order relation on polyhedra Q is geometrical inclusion. 

(b) The least upper bound U is the convex hull. 

(c) The greatest lower bound n is geometrical intersection. 

We will further distinguish the definition of abstract domain operators 
from their actual implementation, which can have bugs. The imple- 
mented version of the operators will be written with a hat: C, U and Fl 
implement the ideal operators 1^, U and Fl, respectively. 

5. inclusion: a polyhedron Pi is included in a polyhedron P2 (noted Pi E -P2) 
if and only if 

Vir,Sat Pi y^Sat P2 y (1) 

Inclusion for constraints Ci = ai ■ x < bi and C2 = 02 • x < 62 is a 
special case which is easy to decide: Ci C C2 holds if and only if there 
exists k > such that k ■ ai = 0,2 and k ■ hi < 62- This latter case is 
thus proven correct directly inside COQ. 

2 Correctness of the Abstract Domain Operators 

Let us now see what needs to be proven for the implementation of each 
operator of an abstract domain so that the correctness of its result can be 
established. 

Inclusion test Pi E P2 ^ Pi E P2 

Convex hull Pi C P1UP2 and P2 Q P1UP2 

Intersection Vx, Sat Pi x A Sat P2 x =^ Sat PinP2 x. 

For now, we will assume a naive implementation of the intersection: 
PinP2 is the union of the constraints of Pi with these of P2, which 
trivially satisfies the desired property M. 

Assignment in a forward analysis, x := e amounts to intersection by the 
equality constraint x' = e (where x' is a fresh variable), projection of x 
and renaming of x' to xjj When analyzing backward, assignment is 
just substitution. 



'^Other polyhedra libraries distinguish invertible assignments (e.g. x := x + 1, more 
generally x' = A- x with A an invertible matrix), which can be handled without projection, 
from non-invertible ones (e.g. x := y + z). Because our library automatically keeps 
a canonical system of equalities, which it uses if possible when projecting, no explicit 
detection of invertibility is needed; it is subsumed by the canonicalization. 



Projection if P2 is the returned polyhedron for the projection of Pi par- 
allel to variables Xj^ , . . . , Xi^ we check that Pi 1^ P2 and that variables 
X jj , . . . , Xjp do not appear in the constraints defining P2 . 

Widening : no correctness check needed. Widening (v) is used to accel- 
erate the convergence of the analysis to a candidate invariant. For 
partial correctness of the analyzer, no property is formally needed of 
the widening operator, since iterations stop when the inclusion test 
reports that an inductive invariant has been obtained. There exist 
formalizations of the widening operator suitable for proving the total 
correctness of the analysis (that is, that it eventually converges to an 
inductive invariant) |15| but we avoided this question by assuming 
some large upper bound on the number of iterations after which the 
analyzer terminates with an error message. 

Remark that we only prove that the returned polyhedron contains the 
polyhedron that it should ideally be (which is all that is needed for proving 
that the results of the analysis are sound), not that it equals it: for instance, 
we prove that the polyhedron returned by the convex hull operator includes 
the convex hull, not that it is the true convex hull. The precision of our 
algorithms (that is, the property that they do not return polyhedra larger 
than needed) is not proved formally; it is however ensured by usual software 
engineering methods (informally correct algorithms, comparing the output 
of our implementation to that of other polyhedra libraries. . . ). 

3 A Posteriori Verification of the Inclusion Test 

We shall now describe a way to ensure the correctness of the inclusion test. 
Recall we represent polyhedra as sets of constraints only. Our certificate 
for proving that a polyhedron P, composed of the constraints Ci, . . . ,Cn 
satisfies a constraint C relies on the following trivial fact: 

Lemma 1 If a point y satisfies a set of constraints {Ci, . . . , Cn}, it satisfies 
any linear positive combination Y17=i ^i^-i with Aj > 0. 

If we can find a constraint C" that is a linear positive combination of Ci , 
. . . ,Cn and such that C QC then it follows that P is included in C. Farkas' 
lemma states that such linear combinations necessarily exist when inclusion 
holds, which justifies our approach. 

The motivation for a posteriori verification of inclusion results stems from 
this formulation: while finding an appropriate linear combination requires 
advanced algorithms, a small program checking that a particular set of Aj's 
entails P^C can easily be proven correct in a proof assistant. We call 
these Aj's the certificate for PQC. 



3.1 A Certificate Checker Certified in Coq 

Our certificate checker has COQ type: 

inclusion.checker {Pi P2 : Polyhedra) {cert : Cert) : Exception (Pi C P2) 

where the type Polyhedra is a simple representation of a polyhedron as 
a list of linear constraints and the type Cert is a representation for in- 
clusion certificates. If a proof of P1QP2 can be built from cert, then 
the inclusion_checker returns it wrapped in the constructor VALUE. How- 
ever cert might be incorrect due to a bug in Q. In this case, the inclusion.checker 
fails to build a proof of Pi CI P2 and returns ERROR. 

When extracting the Ocaml program from the COQ development, proof 
terms are erased and the type of the checker function becomes that which 
would have been expected from a hand-written Ocaml function: [j 

inclusion.checker : Polyhedra — )• Polyhedra — )• Cert — )• bool 

In reality, our implementation is slightly more complicated because the 
untrusted part of our library, for efficiency reasons, operates on fast rational 
and integer arithmetic, while the checker uses standard COQ types that 



explicitly represent integers as a list of bits (see §5.6). 



3.2 A Certificate-Generating Inclusion Test 

Let us now go back to the problem of building a proof of P 1^ C by exhibiting 
an appropriate linear combination. From [4], this can be rephrased as a pure 
satisfiability problem in linear programming: 

(Vy, ^Sat {PnC) y)^ PQC 

This problem can be solved by the simplex algorithm (sl. For this purpose, 
the simplex variant proposed by f9^, designed for SMT-solvers, is particularly 
well-suited. This algorithm only implements the first of the two phases of 
the simplex algorithm: finding a feasible point, that is a point satisfying 
all the constraints of the problem. If there is no such point, a witness of 
unsatisfiability is extracted as a set of mutually exclusive bounds on linear 
terms and suitable Farkas coefficient A, in the same way that blocking clauses 
for theory lemmas are obtained for use in SMT-solving modulo linear rational 
arithmetic. Furthermore, this algorithm is designed for cheap backtracking 
(addition and removal of constraints), which is paramount in SMT-solving 



and also very useful in our application (§5.2) 



^We chose to replace the constructors value and error of the type Exception by 
Ocaml booleans instead of letting the extraction define an Ocaml type "exception" with 
two nuUary constructors due to proof terms being erased. 



Our approach to certificate generation differs from previous suggestions 
[4] where inclusion is first tested by untrusted means, and, if the answer is 
positive, a vector of Farkas coefficients is sought as the solution of a dual 
linear programming problem with optimization, which has a solution, the 
Farkas coefficients, if and only if the primal problem has no solution. Ours 
uses a primal formulation without optimization. 

3.3 From an Unsatisfiability Witness to an Inclusion Certifi- 
cate 

Inclusion certificates are derived from unsatisfiability witness in a way similar 
to [Is]. To illustrate how they are built as part of the inclusion test, a global 
idea of the inner workings of the simplex variant from 19] is needed. We 
insist on the following being a coarse approximation. 

We aim at building, given P non-empty and C, an inclusion certificate 
for PQC, otherwise said P f\C having no solution. P is composed of n 
constraints Ci, . . . , C„ of the form X]?=i ^ij '^j — ^i^ where i is the constraint 

subscript. We refer to C = 60 < J21=i ^Oj ' ^j ^^ ^o- 

Let us start by describing the organization of data. Each constraint Cj is 
split into an equation x'^ = X^^^^ aij ■ Xj and a bound x[ < hi where x[ is a 
fresh variable. For the sake of simplicity, in this presentation, a constraint 
Xi < hi is represented as x[ = Xi /\x'^ <hi ; the actual implementation avoids 
introducing such extra variable. Therefore, each x'^ uniquely identifies Ci by 
construction and the original variables Xi are unbounded. We call basic the 
variables which are defined by an equation (i.e. on the left-hand side, with 
unit coefficient) and non-basic the others. Last, the algorithm maintains a 
candidate feasible point, that is a value for every variable x[ and Xj, initially 
set to 0. 

From this starting point, the algorithm iterates pivoting steps while 
ensuring preservation of the invariant: the candidate feasible point always 
satisfy the equations and the values of the non-basic variables always satisfy 
their bounds (\). At each iteration and prior to pivoting, a basic variable x^ 
is chosen such that its value does not satisfy its bounds. Either there is 
no such x'-, and the candidate feasible point is indeed a solution of P A C, 
thereby disproving P C C; or there is such a basic variable x'^. In this case, 
we shift its value to fit its bounds and we seek a non-basic variable x'^ such 
that its value can be adjusted to compensate the shift: through a pivoting 
step, x[ becomes non-basic, and x^ becomes basic. If there is no such x'„ 
(because all the non-basic variables already have reached their bound), the 
equation which defines x'- exhibits incompatible bounds of the problem and 
is of the form x'^ = YliJ^i -^j " ^'j i^^^Y ^'/s appear in this equation: recall that 
the Xj's are unbounded). We now show how to transform this unsatisfiability 
result into an inclusion certificate. 

Since we supposed that P is not empty, the unsatisfiability necessarily 



involves Cq. Thus, x'q, which represents Co, has a non-zero coefficient Aq in 
the equation. Without loss of generality, we suppose that the incompatible 
bounds involve an upper bound on x'^ and that Aq is positive. The above 
equation can be rewritten so that Xq appears on the left-hand side: 

n 

i=i 

where the lower bound bo < x'q and the upper bound Yll=i K ' ^'j — ^' 
are such that b' < bo. Recall that the x['s were defined as equal to linear 
terms k = X]"=i aij ■ Xj of the constraints Cj. Let us now substitute the x^'s 
by their definition, yielding 

n 

i=i 

Noting that C \slo < bo (since Cq = 60 < ^0 is C), that Yll=i K ' h — ^ ^^^ 
that b' <bo, the A'-'s form an inclusion certificate for PCC. 

4 A Posteriori Verification of the Convex Hull 



We saw in §2 that the result of the convex hull of two polyhedra Pi and P2 
must verify inclusion properties with respect to both Pi and P2. Comput- 
ing P = P1UP2, then Pi !^ P and P2 !^ P and then checking the certificates 
would produce a certified convex hull result, at the expense of two extra in- 
clusion tests. From a development point of view, this is the lightest approach. 
However, careful exploitation of the details of U can save us the extra cost of 
certificate generation, at the expense of some development effort. 

Before delving into the details, let us introduce some more notations for 
the sake of brevity. In this section, a polyhedron P is regarded as a column 
vector of the constraints Ci , . . . , C„ it is composed of. This allows for a 
matrix notation: P = <x \ A- x <b>, where the linear term of Ci is the i 

line of A and the constant of Ci is the i component of b. 

Then, an inclusion certificate, Ai, . . . , A^, for P C C" is a line vector A, 
such that A ■ P = C and C QC' . Now, an inclusion certificate for P^P' 
is a set of inclusion certificates Ai, . . . , A„, one for each constraint C'^ of P'. 
Such a set can be regarded as a matrix F such that 

and P X P C P' 




where the i line of P x P is a constraint C such that C QC'^. We call A a 
Farkas vector and F a Farkas matrix. 



4.1 A Convex Hull Algorithm on Constraints Representa- 
tion 

The convex hull P1UP2 is the smallest polyhedron containing all line segments 
joining Pi to ^2- Thus, a point x of Pi U P2 is the barycenter of a point 
xi in Pi and a point X2 in P2. Exploiting this remark, [3j defined Pi U P2, 
with Pi = {x \ Ai ■ X < bi}, as the set of solutions of the constraints 
^i-a^i < bi A A2-X2 < b2 A X = ai-xi+a2-x'2 A ai+a2 = lA0<aiA0<a2 
using 2n + 2 auxiliary variables aJl, 3J2, ai, 02 where n = |x| is the number 
of variables of the polyhedron. Still following [3j, the variable changes 

02 • X2 remove the non-linearity of the equation 



X^ 



ai ■ xi 



and X2 



X = ai ■ xi + a2 ■ X2. 

The resulting polyhedron can regarded as the 3-block system S,,ar be- 
low. The auxiliary variables x']^ , X2 > «i ; ct2 are then projected out to stick 
to the tuple x of program variables. Therefore, the untrusted convex 
hull operator U mainly consists in a sequence of projections: Pi U P2 = 
proj S^ar (2;i,X2, ai, 02) where 



Jbc 



V 



Aix[ < aibi 



A2X2 < a2b2 



X 


= x\ + 


x'. 


Q] 


+ 02 = 
0< ai 
0<a2 


-- 1 



/ 



4.2 Instrumenting the Projection Algorithm 

Projecting a variable x^ from a polyhedron P represented by constraints 
can be achieved using Fourier-Motzkin elimination (e.g. fsl). This algorithm 
partitions the constraints of P into three sets: E^ contains the constraints 
where the coefficient of Xk is nil, E^ contains those having a strictly positive 
coefficient for x^ and E~ contains those which coefficient for x^ is strictly 
negative. 

Then, the result Pproj of the projection of Xk from P is defined as 



Pproj = proj P Xk = -E9 U ( map elim^^ {Et x E_ 



Xk' 



where E^ x E-. is the set of all possible pairs of inequalities, one element of 

Xk Xk 

each pair belonging to Et and the other belonging to EZ . The elim^^ func- 
tion builds the linear combination with positive coefficients of the members 
of a pair such that Xk has a zero coefficient in the result. 
Illustrating on an example, projecting x from 



P = {y<l,2-x + y<2, 



y < 1} gives 



E^ = {y <1} and E+ X E- = {{2- x + y<2,-x-y< 1)} 
From l-{2-x + y<2) + 2- {-x - y < 1) = -y < 4, we get Pp^oj = {y < 

i,-y<4}. 

Note that every constraint C of -Pproj is either a constraint of P, or 
the result of a hnear combination with non-negative coefficients Ai, A2 of 
two constraints Ci and C2 of P, such that Ai • Ci + A2 • C2 = C. It is 
therefore possible, with some bookkeeping, to build a matrix F such that 
F X P = -Pproj- This extends to the projection of several variables: if 
proj P Xk = Pproj = F X P and proj Pp^j xi = P^^^^ = F' x Pproj, then 
P^^^. = F" X P with F" = F' X F. 

Fourier-Motzkin elimination can generate a lot of redundant constraints, 
which make the representation size of Pproj unwieldy. In the worst case, the 
n constraints split evenly into E^^ and E~^ , and thus, after one elimination, 
one gets n^/4 constraints; this yields an upper bound of n'^'' /4P where p is 
the number of elimination steps. Yet, the number of true faces can only grow 



in single exponential 17, §4.1]; thus most generated constraints are likely to 
be redundant. 

The algorithm inspired from |20| , which we use in practice, adds these 
refinements to Fourier-Motzkin elimination: 

1. Using equalities when available to make substitutions. A substitution 
is no more than a linear combination of two constraints, the coefficients 
of which can be recorded in F. Note that there is no sign restriction 
on the coefficient applied to an equality. 

2. Discarding trivially redundant constraints. The corresponding line F 
can be discarded just as well. 

3. Discarding constraints proved redundant by linear programming, as in 



{5.2 



Note that, since discarding a constraint only adds points to the polyhedron, 
there is no need to prove these refinements to be correct or to provide 
certificates for them. We could thus very easily add new heuristics. 

4.3 On-the-Fly Generation of Inclusion Certificates 

In order to establish the correctness of static analysis, the convex hull 
operator should return a superset of the true convex hull; we thus need 
proofs of Pi E Pi U P2 and P2 !^ Pi U P2 . The converse inclusion is not needed 
for correctness, though we expect that it holds; we will not prove it. A 
certifying operator U must then produce for each constraint C of Pi U P2 a 
certificate Ai (resp. A2) proving the inclusion of Pi (resp. P2) into the single- 
constraint polyhedron C. The method we propose for on-the-fiy generation 
of a correctness certificate is based on the following remark. 



10 



For each constraint C of Pi U P21 the projection operator proj provides a 
vector A such that A x S^ar = C, where S^ar is the system of constraints defined 
An examination of the certificate reveals that A can be spht into three 



!4.1 



m 

parts (Ai, A2, A3) such that Ai refers to the constraints Ai.x[ < aibi derived 
from Pi ; A2 refers to the constraints A2.X2 < 02^2 derived from P2 and A3 
refers to the barycenter part x = x']^ + x 2 A ai + 02 = 1 A < ai A < 02- 
Let us apply the substitution a = [ai/l,a2/0, x'^^/x, X2/O], that characterizes 
the points of Pi as some extreme barycenters, to each terms of the equality 
A X Si,ar = C. This only changes Si,ar'- Indeed, Act = A since A is a constant 
vector and Ca = C since none of the substituted variables appears in C (due 
to projection). We obtain the equality (below) where many constraints of 
SbarO' became trivial. 



(Ai,A2,A3) X 



/ 


^ix < 6] 


\ 












0< 










\ 


X = X 

1 = 1 

0< 1 
0< 


/ 



c 



This equality can be simplified into Ai x {Aix < 6i) + A(0 < 1) = C 
where A is the third coefficient of A3. This shows that Ai is a certificatq^ 
for Pi C C. The same reasoning with a = [ai/0, a2/l, af'i/O, X2/X] shows that 
A2 is a certificate for P2 C C. 

5 Notes on the Implementation 

The practical efficiency of the abstract domain operators is highly sensitive 
to implementation details. Let us thus describe our main design choices. 

5.1 Extending to Equalities and Strict Inequalities 

Everj^hing we discussed so far deals with non-strict inequalities only. The 
inclusion test algorithm however complements such non-strict inequalities, 
which yields strict ones. Adaptation could have been restricted to the simplex 
algorithm on which the inclusion test relies, and such an enhancement is 
described in [9l]. We have however elected to add full support for strict 
inequalities to our implementation. Once the addition of two constraints 
has been defined, almost no further change to the algorithms we discussed 
previously was needed. 

Proper support and use of equalities was more involving. As [20] points 
out, equalities can be used for projecting variables. Such substitutions do not 



^The shift A of the bound is lost and wiU be computed again by our COQ-certified 
checker. 



11 



increase the number of constraints, contrary to Fourier-Motzkin elimination. 
We ended up splitting the constraint set into a set of equalities, each serving 
as the definition of a variable, and a set of inequalities in which these 



variables have been substituted by their definitions. Minimization (see §5.2) 
was augmented to look for implicit equalities in the set of inequalities. Last, 
testing inclusion of P in C was split into two phases: substituting in C the 
variables defined by the equalities of P and then using the simplex-based 
method described earlier without putting the equalities of P in, which reduces 
the problem size. 

Inclusion certificates were adapted for equalities. If P C C, with C = a- 
X = b, cannot be proven using a linear combination of equalities, it is split 
as {a • X < b, a ■ X > b} and P is proven to be included in each separately. 

5.2 Minimization 



The intersection PinP2 is a very simple operation. As §2 described, a naive 
implementation amounts to list concatenation. However, some constraints 
of Pi may be redundant with constraints of P2- Keeping redundant con- 
straints leads to a quick growth of the representation sizes and thus of 
computation costs. In addition, one condition for the good operations of 
widening operators on polyhedra is that there should be no implicit equality 
in the system of inequalities and no redundant constraint [2] . 

It is therefore necessary to minimize the size of the representation of 
polyhedra, that is, removing all redundant constraints, and to have a system 
of equality constraints that exactly defines the affine span of the polyhedron. 
We call Pmin the result of the minimization on P. The correctness of the 
result is preserved as long as Pmin is an over-approximation of P, which 
means P I^Pmin- 

First, we check whether P has points in it using the simplex algorithm 



from §3.3 If P is empty, _L is returned as the minimal representation. The 
certificate is built from the witness of contradictory bounds returned by 
the simplex algorithm. It is a linear combination which result is a trivially 
contradictory constraint involving only constants (e.g. < —1) and which, 
in other words, has no solution. 

The next step is implicit equality detection. It builds on a-x < 6 A a-x > 
b ^ a ■ X = b. For every C- = a ■ x <b oi P (by definition P C C-), we 
test whether P C C- = a- x > b. If the inclusion holds, the certificate of the 
resulting equality is composed of a linear combination yielding C- and a 
trivial one, 1 • C- , yielding C- . Once this is done, the representation of P 
can be split into a system of equalities Pe and a system of inequalities Pi 
with no implicit equality. Pe is transformed to be in echelon form using 
Gaussian elimination, which has two benefits. First, redundant equations 
are detected and removed. Second, each equation can now serve as the 
definition of one variable. The so-defined variables are then substituted in Pj, 
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yielding P^. Although our implementation tracks evidence of the correctness 
of this process, it should be noted that the uses of equalities decribed above 
are standard practice. 

At this point, if redundancy remains, it is to be found in P^ only. It 
is detected using inclusion tests: for every C £ Pi, if P/ \ {C} E-P/, C is 
removed. Removing a constraint is, at worst, an over-approximation for 
which no justification needs to be provided. 

All that we describe above involve many runs of the simplex algorithm. 
The key point which makes this viable in practice is the following: they are 
all strongly related and many pivoting steps are shared among the different 



queries. We described (§3.3) the data representation used by the simplex 
variant we use: it splits each constraint of P in linear term and bound 
by inserting new variables. These variables can have both an upper and 
a lower bound. Let us now illustrate the three steps of minimization on 
constraint C = a ■ x < b, split as x' = a ■ x and 3? < h. The first step, 
satisfiability, solves this very problem. Then, implicit equalities detection 
checks whether x' = a ■ x and x' < 6 is unsatisfiable. Last, redundancy 
elimination operates on x' = a ■ x and x' > b. 

For all these problems, we only changed the bound on x' , without ever 
touching either the constraint x' = a- x or the other constraints of P. These 
changes can be done dynamically, while preserving the simplex invariant 



(I of §3.3), by making sure that the affected x' is a basic variable. This 
remark, once generalized to a whole polyhedron, enables the factorization of 
the construction of the simplex problem. Actually, it is only done once for 
each minimization. It is also hoped that the feasible point of one problem is 
close enough to that of the next problem, so that convergence is quick. 
Minimization also plays an important role in the convex hull algorithm. 



We mentioned (§4.2) that projection increases the representation size of 
polyhedra and described some simple counter-measures from [26] . When 
projecting a lot of variables, as is done for computing the convex hull of two 
polyhedra, each redundant constraint can trigger a lot of extra computation. 
Applying a complete minimization after the projection of each variable 
mitigates this. More precisely, only the third of the steps described above is 
used: projection cannot make a non-empty polyhedron empty and it cannot 
reduce the dimension of a polyhedron, no implicit equality can be created. 

5.3 A More Detailed Intuition on Bookkeeping 



We mentioned in §3.2 and §4 that simple bookkeeping makes it possible to 
build inclusion certificates. We now give a more precise insight on what is 
involved, on the example of the projection. 

The main change is an extension of the notion of constraint, which is 
now a pair (/, C) of a certificate fragment and a linear constraint as we 
presented them so far. A certificate fragment / is a list of pairs (n^, idi), ni 
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being a rational coefficient and idi a natural number uniquely identifying 
one constraint of P. The meaning of / is the following 

y^ rii ■ Cidi = c, with Cid^ £ P and (nj, idi) £ f 



The elirrig^, function introduced in|<ii4.2|is extended to take two extended 



constraints (/i,Ci) and (/2,C2), and return an extended constraint (/, C). 
Recall that the original elim^^ chooses Ai and A2 such that the coefficient 
of Xk in the resulting C is nil. The extended version returns (Ai • /i @ A2 • 
/2, Ai • Ci + A2 • C2), where @ is the list concatenation operator and Aj • fi is 
a notation for: 

map {fun {n, id) — )• (Aj • n, id)) fi 

The certificate fragment keeps track of how a constraint was generated 
from an initial set of constraints. For a single projection proj P x^, the 
fragments are initialized as [(1, idc)] for every constraint C before the actual 
projection starts. For a series of projection as done for the convex hull, the 
initialization takes place before the first projection. 

5.4 Polyhedron Representation Invariants 

The data representation our implementation uses for polyhedra satisfies a 
number of invariants which relate to minimality. 

(1) There is no implicit equality among the inequalities. 

(2) There is no redundant constraint, equality or inequality. 

(3) In a given constraint, factors common to all the coefficients of variables 
are removed. 

(4) Each equality provides a definition for one variable, which is then substi- 
tuted in the inequalities. 

(5) Empty polyhedra are explicitly labeled as such. 



I (3) I helps keeping numbers small, hopefully fitting machine representation, 
resulting in cheaper arithmetic. |(1)| implies in particular that if an implicit 
equality is created when adding a constraint C to a polyhedron P, then C is 
necessarily involved in that equality. It follows that the search for implicit 
equalities can be restricted to those involving newly added constraints. 



Because of (2) the same holds for redundancy elimination: if C is shown to 
be redundant, P remains unaffected by the intersection. Furthermore, [(4)] 
allows for the reduction of the problem dimension when testing for Pi C ^2- 
Once the same variables are substituted in Pi and P2, only the inequalities 
need to be inserted in the simplex problem. Last, |(3)| and |(4)| give a canonical 



14 



form to constraints, which make syntactic criteria for deciding inclusion 
of constraints more powerful. These criteria, suggested by [20j, are used 
whenever possible in the inclusion test and the projection. 

5.5 Data Structures 
5.5.1 Radix Trees 

Capturing linear relations between program variables with polyhedra gener- 



ally leads to sparse systems, as noted by 20 . Our implementation uses a tree 
representation of vectorq^ where the path from the root to a node identifies 
the variable whose coefficient is stored at that node. This offers a middle 
ground between dense representation, as used by other widely- used imple- 
mentation of the abstraction domain of polyhedra, and sparse representation 
which makes random access costly as sparsity diminishes. 

5.5.2 Numbers Representation 

Rational vector coefficients can grow so as to overflow native integer repre- 
sentation during an analysis. Working around this shortcoming requires the 
use of an arithmetic library for arbitrarily large numbers. This has a serious 
impact on overall performance. Our implementation uses the ZArith[14] 
OCAML front-end to GMPf22]. ZArith tries to lower the cost of using GMP 
by using native integers as long as they don't overflow. 

Our experiments show that, in many practical cases, extended precision 
arithmetic is not used. This echoes similar findings in SMT-solvers such as 



Z3 or OpenSMT 19 : in most cases, extended precision is not used, thus the 
great importance of an arithmetic library that operates on machine words 
as much as possible, without allocating extended precision numbers. In the 
case of polyhedra, however, the situation occasionally degenerates when the 
convex hull operator generates large coefficients. 

The extracted Ocaml code of inclusion.checker does not use this 
efficient representation. Because of the need for correctness of computations, 
the checker instead uses the COQ representation of numbers (lists of bits), 
which is inefficient on numerical computations. Alternatively, assuming trust 
in ZArith and GMP, it is possible to configure the COQ extractor to base 
the checker on ZArith. 

5.6 A Posteriori Certification vs. Full COQ-Certified Devel- 
opment 

Even though our library is planned to be used in a COQ-certified analyzer, 
we preferred a posteriori certification over a fully COQ-certified development. 



*The idea was borrowed from |4|. 
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Keeping COQ only for the development of checkers of external computa- 
tions reduces the development cost and reconciles efficiency of the tool and 
confidence in its implementation through certificates. 

First, it reduces the proof effort: verifying that a guess is the solution to a 
problem involves weaker mathematic arguments than proving correctness and 
termination of the solver. To illustrate the simplicity of our COQ development, 
Figure [T] shows some excerpts which are self-explanatory. The last function, 
inclusion.checker, is representative of the difficulty of the proofs. This 
function is close to its extraction in Ocaml except that it returns either an 
ERROR or a proof of Pi C P2 wrapped in the value constructor (Line 38). In 
the case where Pi is an empty polyhedron (established by eproof) the proof 
of inclusion in P2 is built from that proof of emptiness. The missing proof 
of Line 38 is done in the interactive prover (Lines 43-45) and automatically 
placed in the function. It consists in an induction on the list of constraints of 
P2 that shows that the empty polyhedron Pi is included in every constraint 
ofPs. 

Our external library acts as an oracle: it efficiently performs the operation 
and returns a certificate which serves two purposes: it can be used to check 
the correctness of the computations but it is also a short cut toward the 
result. For instance, the convex hull Pi U P2 is easy to obtain from the 
complete inclusion certificates (Pi,Ai) related to Pi o r (F 2, X2) related to 



P2. Indeed, Pi U P2 = Pi • Pi + Ai = P2 • P2 + A2 (see ^4.31 ). This way, the 
expensive computations that involve numerous calls to our simplex algorithm 
are done by our Ocaml implementation using ZArith and the result is 
refiected in COQ at the cost of just a matrix product using the COQ-certified 
representation of numbers. If we work in such a manner, we never actually 
have to transfer polyhedra from the untrusted to the trusted side. 

From a general point of view, splitting a tool into an untrusted solver and 
a correctness checker makes it more amenable to extensions and optimizations. 
A posteriori certification has a cost each time the correctness of a result 
needs to be proved (only during the last phase of the analysis to ascertain the 
stability of the inferred properties). However, it allows optimizations whose 
correctness would be difficult to prove and usage of untrusted components 
(e.g. GMP). 

6 Experimental Results 

In order to evaluate the viability of our solution, we compared experimentally 
our library (referred to as Libpoly) with mature implementations. 

In addition to the efficiency of the polyhedra computation, we wished 
to measure the cost of the inclusion checker. Our approach guarantees that, 
if our certificate checker terminates successfully on a given verification, the 
result of the operation which produced the certificate is correct. However, 
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From module LinearCsrt: 

Record LinearCstr: Set := mk {coefs: Vec; cmp_op: Cmp; bound: Num} . 

Definition Sat (c: LinearCstr) (a;: Vec) : Prop : = 

denote (Vec.eval (coefs c) ;c) (cmp_op c) (bound c) . 

From module List: 

Inductive Forall (A : Type) (pred : A — > Prop) : list A — > Prop : = 

I Forall_nil: Forall pred nil 

I F0RALL_C0NS: V (.X : k) (l:list A), 

pred X — >■ Forall pred 1 — >■ Forall pred (a; : : 1) 

From module Polyhedra: 

Definition Polyhedra : Set := list (id * LinearCstr). 

Definition Sat (P:Polyhedra) (xr.Vec) : Prop : = 

List. Forall (fun c => LinearCstr .Sat (snd c) x) P. 

Definition Incl (P: Polyhedra) (C: LinearCstr) : Prop : = 

V 2;:Vec, Sat P a: — > LinearCstr .Sat C x. 

Definition (infix C) (Pi P2 : Polyhedra) : Prop := 

V x:Vec, Sat Pi a; ^ Sat P2 x. 

Definition CertOneConstraint : Set := list (id * Num) 

Inductive Cert : Set := 
I incl: list (id * CertOneConstraint) -> Cert 
I empty: CertOneConstraint -> Cert. 

Lemma Empty_is_included: V (P:Polyhedra) (C : LinearCstr) , 
(Empty P) -^ (Incl P C) . 

Definition inclusion_checker (Pi P2 : Polyhedra) (cert:Cert) : Exc(PiCP2). 

refine ( match cert with 

I INCL icert => checklnclusion Pi P2 icert 

I EMPTY ecert => match (checkEmptyness Pi ecert) with 

I VALUE eproof => value _ ^— missing proof 

I ERROR => error 

end 
end 
) . The missing proof is provided by the following proof script: 
induction P2 with IH; 
exact (List .F0RALL_NIL _ _) ; 

exact (List .F0RALL_C0NS _ _ c _ (Empty_is_included Pi (snd c) eproof) IH) , 
Defined. 



Figure 1: Excerpts of our COQ-certified inclusion checker 
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this assertion currently only applies to the polyhedra as known to the COQ 
checker: a translation occurs between the Ocaml representation of numbers, 
ZArith, and their representation in the COQ language as lists of bits. This 
means that the checker has to compute on this inefficient representation, and 
thus we wished to ascertain whether the cost was tolerable^ 

The best approach to evaluating Libpoly would have been to rely on 
it for building a complete static analyzer. Although this is our long-term 
goal, a less demanding method was needed for a more immediate evaluation. 
We chose to compare computation results from LiBPOLY to those of widely 
used existing implementations of the abstract domain of polyhedra: the 
NewPolka library and the PPL. More precisely, we used them through 
their Apron front end [II] . 

6.1 The Method 



As [16] points out, randomly-generated polyhedra do not give a faithful 
evaluation: a more realistic approach was needed. Because of the lack of a 
static analyzer supporting both Apron and Libpoly, we carried out the 
comparison by logging and then replaying with Libpoly the abstract domain 
operations done by the existing Pagai analyzer |10| using Apron. 

Technically, logging consists in intercepting calls to the Apron shared 
library (using the wrap functionality of the GNU linker Id), analyzing the 
data structures passed as operands and generating equivalent Ocaml code 
for Libpoly. NewPolka and PPL results are logged too, for comparison 
purposes. At the end of the analysis, the generated Ocaml code forms a 
complete program which replays all the abstract domain operations executed 
by the NewPolka library or the PPL on request of the analyzer. 

The comparison was done for the following operations: parallel assign- 
ment, convex hull, inclusion test and intersection on the analysis of the 
following programs: 

1. bf : the Blowfish cryptographic cipher 

2. bz2: the bzip2 compression algorithm 

3. dbz2: the bzip2 decompress algorithm 

4. jpg: an implementation of the jpeg codec 



""An alternative would be to map, at checker extraction time, COQ numbers to ZArith 
numbers, at the expense of having both ZArith and GMP in the trusted computing 
base. One may consider that we already make assumptions about ZArith and GMP: we 
assume they respect memory safety, and thus will not corrupt the data of the Ocaml code 
extracted from COQ, or at least that, if they corrupt memory, they will cause a crash in 
the analyzer (probably in the garbage collector) instead of a silent execution with incorrect 
data. This seems a much less bold assumption than considering that they always compute 
correctly, including in all corner cases. 



18 



5. re: the regular expression engine of GNU awk 

6. f oo: a hand-crafted program leading to polyhedra with many con- 
straints, large coefficients and few equalities 

6.2 Precision and Representation Size Comparison 

The result of each operator we evaluated is a well-defined geometrical object. 
For every logged call, the results from NewPolka, PPL and Libpoly were 
checked for equality (double inclusion). The certificates generated by Libpoly 
were then systematically checked. Furthermore, polyhedra have a minimal 
constraints representation, up to the variable choices in the substitutions of 
equalities. It was systematically checked whether Libpoly, NewPolka and 
the PPL computed the same number of equalities and inequalities. In all 
the cases we tried, the tests of correctness and precision passed. It is to be 
noted that the PPL does not systematically minimize representations: its 
results often have redundant constraints|£] 

Besides giving confidence in the results computed by Libpoly, ensuring 
that our results are identical to those of NewPolka or the PPL lead us to 
believe that the analyzer behavior would not have been very different, had it 
used the results from Libpoly. There is no noticeable difference between 
the analyses carried out using NewPolka and the PPL. 

6.3 Timing Measurements 

Timing measurements were made difficult because of the importance of the 
state of polyhedra in the double representation NewPolka and the PPL 
use. We were concerned that logging and replaying as described above 
would be unfair towards these libraries, since it would force the systematic 
recomputation of generator representations that, in a real analyzer, would 
be kept internally. We thus opted for a different approach. 

We measured the timings for NewPolka and the PPL directly inside 
Pagai by wrapping the function calls between calls to a high precision timer. 
We made sure that the overhead of the timer system calls was sufficiently 
small so as to produce meaningful results. For Libpoly, timing measurements 
were done during the replay and exclude the time needed to parse and rebuild 
the operand polyhedra. 

We present two views of the same timing measurements, carried out 
on the programs introduced in §6.1 Table [T] gives, for each benchmark 



program, the total time spent in each operation of the abstract domain. 
Such a table does not inform us of the typical distribution of problem sizes 



^This is due to the lazy-by-default implementation of the operators of the PPL. Since 
support for the eager version of the operators has been deprecated in and is being removed 
from the PPL (see [2^, § A Note on the Implementation of the Operators), we could not 
configure the library to have the same behavior as NewPolka. 
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Table 1: Timing comparison between NewPolka (N), PPL (P), Libpoly 
(L) and Libpoly with certificate checker (C): total time (in milliseconds) 
spent in each of the operations; trivial problems are excluded. 



prog. 



assignment 
N P L 



N 



convex hull 
P L 



C 



inclusion 
N P L C 



intersection 
N P L 



bf 

bz2 

dbz2 

JPg 
re 

f 00 



3.7 11.4 0.5 

14.6 54.1 2.9 
1618 4182 83.8 

23.7 68.3 3.8 
5.7 17.2 0.7 
9.2 14.8 8.5 



3.2 1.2 2.7 2.8 

23.5 11.5 66.8 68.7 

1393 231.9 532.8 535.3 

28.2 7.5 24.0 24.9 

20.2 8.4 17.9 19.2 

4.2 0.6 941.8 943.7 



0.2 0.4 0.10.1 
1.6 2.8 0.7 1.2 
32.3 35.6 2.1 3.6 
1.2 1.8 0.5 0.8 
1.1 1.3 0.5 0.7 
0.2 0.2 0.9 0.9 



10.7 13.4 
52.3 61.1 
1687 1815 
39.7 51.0 
37.3 47.2 
6.7 7.1 



1.2 
7.9 
28.3 
6.0 
3.3 
5.5 



and the relationship between problem size and computation time, thus we 
compiled Table [2] which shows computation times aggregated according to 
the "problem size" , defined as the sum of the number of constraints of all 
the operands of a given operation. 

For the assignment and the convex hull, all the constraints of the two 
operands are put together after renaming and many projections follow. The 
inclusion test Pi ^P2, in the worst case, solves as many linear programming 
problems as there are constraints in P21 but each is of size the number of 
constraints of Pi + 1. Last, the intersection operator minimizes the result of 
the union of the sets of constraints. Note that the sums in Table [T] exclude 
operations on trivial problems of size zero or one. 

The presented results show that Libpoly is efficient on small problems. 
Yet, the performance gap between Libpoly and the other implementations 
closes on bigger problems. This is especially true for the convex hull, which 
is a costly operation in the constraint representation. At least part of the 
difference in efficiency on small problems can be explained by the generality 
Apron provides: it provides a unified interface to several abstract domains 
at the expense of an extra abstraction layer which introduces a significant 
overhead on small problems. 

More generally, the use of ZArith in Libpoly is likely to lower the cost 
of arithmetic when compared to NewPolka and the PPL, which use GMP 
directly. The f 00 program illustrates this: the analysis creates constraints 
with big coefficients, likely to overflow native number representation. However, 
precise measurement of the effect of using ZArith would be a hard task. 

Last, Table [T] seems to show that problems are most often of rather small 
size, but this may well be due to our limited experimentation means. 

In spite of the shortcomings of our evaluation method, these results seem 
promising for a constraints-only implementation of the abstract domain of 
polyhedra. Some progress still needs to be made on the convex hull side (see 
§7). It is also interesting to notice the performance differences between the 
NewPolka and the PPL, despite their design similarities; we ignore their 
cause. 
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Table 2: Timing comparison between NewPolka (N), PPL (P) and Libpoly 
(L). Computation times (in milliseconds) are aggregated according to opera- 
tion and problem size, (n) is the total number of problems of the size range 
in the benckmarks. 



problem size 


0-1 


2-5 


6-10 


11-15 


16-20 


21-25 


26-30 


31+ 


N 


33.8 


601.8 


385.4 


20.9 


78.3 


537.4 


59.5 


13.1 


assignment P 


47.5 


1176 


519.7 


87.4 


247.6 


2111 


81.7 


77.9 


L 


1.1 


6.6 


14.3 


10.7 


5.2 


39.2 


15.2 


11.6 


n 


539 


667 


381 


58 


64 


480 


30 


16 


N 


687.9 


679.7 


434.1 


119.5 


68.8 


37.9 


6.4 


3.5 


convex hull P 


167.5 


141.0 


68.4 


22.8 


16.8 


9.2 


1.9 


0.9 


L 


7.0 


57.1 


133.7 


131.2 


1050 


106.4 


50.1 


27.8 


n 


3354 


3373 


1092 


354 


135 


65 


U 


7 


N 


7.2 


9.7 


9.7 


3.3 


5.8 


4.0 


4.0 





inclusion P 


6.5 


12.8 


10.6 


4.2 


7.0 


3.9 


3.4 





L 


0.6 


1.6 


1.3 


0.5 


1.0 


0.3 


0.1 





n 


U82 


1881 


673 


277 


111 


52 


17 


4 


N 


1389 


1752 


52.3 


27.4 


1.3 








intersection P 


1933 


1740 


158.6 


91.4 


4.8 








L 


35.0 


30.9 


18.4 


8.8 


0.6 








n 


11458 


4094 


322 


156 


6 












6.4 Certificate Cliecking Overhead 

The certificate checking overhead shown in Table [T] includes the translation 
between Ocaml and COQ representations. Inside a certified static analyzer, 
this overhead could be reduced by only transferring the certificates, as opposed 
to the full polyhedra, and using them to simulate the polyhedra computations, 
without bothering to check after every call that the polyhedron inside the 
Ocaml library corresponds to the one inside the certified checker. In addition 
to translation costs, there is the general inefficiency of computations on COQ 
integers, which are represented as lists of bits; this is considerably more 
expensive than using native integers, or even arrays of native integers as 
GMP would do. 

However, it should be noted that the checking of inclusion certificates 
occurs only during the final step of the certified static analysis which consists 
in verifying that the inferred invariant candidates are indeed inductive 
invariants for the program. 

Last, the overhead of certificate checking is relatively greater for inclusion 
than for convex hull. Although the actual checking burden is bigger for the 
convex hull, due to certificate composition densifying the resulting certificate, 
the inclusion test algorithm is much cheaper than the convex hull in terms of 
computations. More precisely, the convex hull algorithm involves inclusion 
tests as part of representation minimization. 
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7 Conclusions 

The previous sections demonstrated that a reahstic implementation of the 
abstract domain of polyhedra can be certified using a posteriori verification 
of results. This approach has a key benefit: the time-consuming develop- 
ment inside the COQ proof assistant is reduced to the bare minimum. A 
tight integration of the certification concern enables on-the-fly certification 
generation as a by-production of the actual computations, thereby making 
the associated cost negligible. The same procedures can be used for fixed 
point iterations (with certificate generation turned off for efficiency) and for 
fixed point verification (with certificates generated and checked). 

The complete implementation which has been developed operates only on 
a constraints representation of polyhedra; our motivations for this choice were 
the ease of generation of certificates as well as the absence of combinatorial 
explosion on common cases such as hypercubes. This is made possible through 
careful choice of data structures and exploitation of recent algorithmic 
refinements |20[ M . Possible future developments include designing efficient 
techniques for generating Farkas certificates for a library based on the double 
representation (generators and constraints) and providing heuristics for 
choosing when to operate over constraints only and when to use the double 
representation. 

Prior to this, however, there remains room for both enhancement and 
extension of our current implementation. A simple enhancement would be to 
have both an upper and a lower bound for linear terms, which would further 
condense the representation of polyhedra. The implicit equality detection 
algorithm could be made less naive by exploiting the fact that a point in a 
polyhedron P which has implicit equalities Ei necessarily reaches the bounds 
of the inequalities involved in the proof of PQEi. 

Finally, our library is planned to be part of a certified static analyzer, 
such as the one being built in the Verasco project. Beyond a certified 
implementation of the abstract domain of polyhedra, our library could also 
serve to verify the numerical invariants discovered by untrusted analysis 
using a combination of abstract domains (intervals, octagons, ... which are 
special cases of polyhedra) . The discovered invariants could be stored in the 
form of polyhedra and the verification of their stability could be done with 
our certified library. Currently, our polyhedron library only deals with linear 
constraints, but a general-purpose analyzer has to handle nonlinearity. Our 
library should therefore include linearization techniques p3j at the condition 
that these be proven correct. 

Acknowledgements We would like to thank Bertrand Jeannet for his 
advice on proper ways to evaluate Libpoly against his NewPolka library. 
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